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. ■ Abstract 
<_ 

| 7j . Mixed dictionaries generated by cosine and B-spline functions are considered. It is shown that, 

(-^ I by highly nonhnear approaches such as Orthogonal Matching Pursuit, the discrete version of 

I the proposed dictionaries yields a significant gain in the sparsity of an image representation. 

1 Introduction 

> , 

O ! Sparse representation of information is a central aim of data processing techniques. An usual 
first step of image processing applications, for instance, is to map the image onto a transformed 
space allowing for the reduction of the number of data points representing the image. Currently 
Q> ■ the most broadly used transforms for performing that task are the Discrete Cosine Transform 
^ . (DCT) and Discrete Wavelet Transforms (DWT). An important reason for the popularity of 
Q ! both these transforms is the viability of their fast implementation. However, since parallel 
processing is becoming more powerful and accessible, alternative approaches for signal repre- 
sentation are being given increasing consideration. Emerging techniques address the matter in 
the following way: Given a signal / G find the decomposition / = Xli^li '^i'^^i 5 where vectors 
■ v^. G M^, i = 1, . . . , M, usually called atoms, are a subset of a redundant set called a dictionary. 
Approximations of this type are highly nonlinear and are said to yield a sparse representation 
of the signal / in terms of M atoms if M is considerably smaller than N. Available method- 
ologies for nonlinear approximations are known as Pursuit Strategies. This comprises Bases 
Pursuit [5] and Matching Pursuit like greedy algorithms [9], including Orthogonal Matching 
Pursuit (OMP) and variations of it [3,4,11,12]. Another concern inherent to highly non lin- 
ear approximations is the design of suitable dictionaries for representing some class of signals. 
Dictionaries arising by merging orthogonal bases are theoretically studied in [6,7]. From a 
different perspective, approaches for learning dictionaries from large data sets are considered 
in [8, 10]. In this communication we present an alternative construction of dictionaries for rep- 
resenting natural images. The proposed dictionaries are a mixture of discrete cosine and spline 
based dictionaries. We have found by a good number of examples, some of which are presented 
here, that the resulting dictionary renders a considerable gain in sparsity in comparison to fast 
transforms such as DCT and DWT, at acceptable visual level (PSRN 40 dB). 

The paper is organized as follows: In Sec. [2] we introduce the discrete B-spline based dic- 
tionaries which together with the discrete cosine ones form the large mixed dictionary we are 



1 




Figure 1: Splines taken from dictionaries spanning the same space. Linear splines (left) Cubic 
splines (right) 

proposing. In Sec. [3] we discuss the implementation of the OMP approach that we have used. 
The details of the actual process for dealing with images are given in Sec. H] where results illus- 
trating the capability of the proposed dictionaries to yield sparse representations by nonlinear 
approaches are provided. 

2 B-spline based dictionaries 

The discrete dictionaries we discuss here are inspired by a general result holding for continuous 
spline spaces. Namely, that spline spaces on a closed interval can be spanned by dictionaries of 
B-splines of broader support than the corresponding B-spline basis functions [2]. 

A partition of an interval [c, d] is a finite set of points A := {xi}^J^^, N G N such that c = 
Xo < a;i < ■ ■ ■ < xat < x^+i = d, which generates N subintervals = [xi, Xj+i), i = 0, . . . , N — 1 
and In = [s^at, a^A^+i]- Representing by 11^ the space of polynomials of degree smaller than or 
equal to m G No = N U {0} and as C" the space of functions having continuous derivatives 
up to order m (with the space of continuous functions) the spline space of order m > 2 on 
[c, d] , with single knots at the partition points, is define as 

SmiA) = {feC"'-%d] : /I,, Gn^„i,z = 0,...,iV}, 

where /|/. indicates the restriction of the function / to the interval Jj. In the case of equally 
spaced knots the corresponding splines are called cardinal. Moreover all the cardinal B-splines 
of order m can be obtained from one cardinal B-spline B{x) associated with the uniform simple 
knot sequence 0,1, ... ,m. Such a function is given as 

B^x) = - J2{-ly^]{x-^)^-\ (1) 
1=0 ^ ^ 

where {x — i)+~^ is equal to {x — i)'^~^ ii x — i > and otherwise. We shall focus on the 
particular cases corresponding to m = 2 and m = 4. For m = 2 the cardinal spline space 
5*2 (A) is the space of piece wise linear functions and can be spanned by a linear B-spline basis 
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Figure 2: Tree splines taken from dictionaries spanning the same space. Linear splines (left) 
Cubic splines (right) 

arising by translating the prototype function known as 'hat' function. The first 3 functions in 
the left graph of Fig. [T] are 3 consecutive linear B-spline basis functions. The 3 middle functions 
in the same graph are linear B-spline functions of broader support taking from a dictionary 
spanning the same space as the basis. The last 3 functions are taking from another dictionary 
for the same space. Details on how to build these dictionaries are given in [2]. The basis and 
dictionary functions equivalent to the ones in the left graph of Fig.[Tl but for cubic spline spaces 
corresponding to m = 4, are given in the right graph of the same Figure. 

For constructing redundant dictionaries suitable for processing images by nonlinear tech- 
niques we need to make sure that the dictionaries can be processed with digital computers 
having the existing memory capacity. Thus, we need to a)discretise the functions to obtain ad- 
equate Euclidean vectors and b)restrict the functions to small intervals allowing for processing 
the images in small blocks. We carry out the discretization by taking the value of a prototype 
function only at the knots (cf. small circles in graphs Fig. [2]) and translating that prototype 
one sampling point at each translation step. In regard to the boundaries one may take different 
routes: A possibility is to adopt periodicity (cyclic boundary conditions) and other apply the 
'cut off' approach and keep all the vectors whose support has nonzero intersection with the in- 
terval being considered. The former would leave a basis for the corresponding Euclidean space 
and the later a redundant dictionary. 

Remark 1. We notice that by the proposed discretization the hat B-spline basis for the cor- 
responding interval becomes the standard Euclidean basis for either boundary conditions. By 
discretizing the hats of broader support the samples preserve the hat shape. 

Obviously for a finite dimension Euclidean space we can construct arbitrary dictionaries. 
In particular, different B-spline based dictionaries each of which comprising vectors of different 
support. Furthermore, we can include vectors of different support by merging dictionaries. 
There is, of course, a compromise between redundancy and complexity that needs to be consid- 
ered. The discussion of such a tradeoff is postponed to Sec. HI where the numerical examples 
are described. 

From discrete unidimensional B-spline based dictionaries we obtain bidimensional ones sim- 
ply by taking tensor product. Actually in Sec. IHwe consider a dictionary consisting of unidi- 
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mensional cosine and B-spline based vectors, of redundancy approximately five, and build the 
bidimensional one by taking the tensor product of the whole dictionary with itself. 



3 Implementation of the greedy algorithm OMP 

The OMP technique [11] is an adaptive greedy strategy for selecting atoms which evolves as 
follows: Let / G be a given signal and {f a given redundant dictionary. Setting = f 
at iteration k + 1 the OMP algorithm selects the atom, f j^.^^ say, as the one minimizing the 
absolute values of the inner products {vi, R'^), i = 1, . . . , L, i.e., 

k 

= arg max | {vi, R^) \ , where R'' ^ f -V] c^v^. , (2) 

1=1 

and J is the set of indices labehng the dictionary's atoms. The coefficients , i = 1, . . . , A; in the 
above decomposition are such that ||/ — is minimum, which is equivalent to requesting 

R^ = Pvkfi where Py^ is the orthogonal projection operator onto Vk = spsin{v£^}f^i. Wc 
base our implementation for determining the coefficients q, i = l,...,k on Gram Schmidt 
orthogonalization with re-orthogonalization, and recursive biorthogonalization. Basically, at 
each iteration we update the vectors 

where v^^l = with g^+i = Vi^^^^ - PvkQk+i and qi = v^^. One reorthogonahzation 

step implies to recalculate Qk+i as Qk+i = Qk+i — Pv^lk+i- The projector Py^. is here computed 
as Pyfc = QkQk where the fc-columns of matrix Qk arc the vectors (/i/||q't!|, i = 1, .... A; and 
Ql indicates the transpose conjugate of Qk- However, to calculate the coefficients of the linear 
superposition we express the projectors as Py^ = ^kBl where the A;-columns of matrix are 
the selected vectors and the /c-columns of matrix Bk are the vectors vf, i — 1, . . . , /c. Thus, 
the required coefficients arise from the inner products = {vf, f), i = 1, . . . ,k. Details on 
this type of implementation are given in [4, 12] and the code can be found at [1]. Moreover, as 
will be discussed in the next section, the fact that we deal with dictionaries involving DC and 
supported atoms reduces the general complexity of the OMP method. 

4 Sparse image representation by Discrete Cosine and 
B-spline based dictionaries 

Here we present the examples of the gain in sparsity that are achieved by using dictionaries 
which are the union of Discrete Cosine an B-Sphne based dictionaries. In order to apply the 

OMP approach on an image using dictionaries of theses types is necessary to divide the image 
into blocks. The size of all the test images we consider is 512 x 512 and wc divide each image 
into blocks of 16 x 16 pixels. For approximating each block we first construct the dictionaries 
I>i and I>i, = 2, . . . , 4 defined as follows 

• Discrete Cosine Dictionary 



P.^{c.co.( "ft'-^''-'> U^l,...,iV}U 



I M 
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with Cj = ^i = l,...,M and Ci = —k^ a normalization factor. 
Discrete B-Spline based dictionaries 



P, = {^.i?^(j-z)|L;j = l,...,L}fl^i, 

where the notation Bm{j — i)\L indicates the restriction to be an array of size L, indices 
k = 2, ... ,4 label the dictionaries of different support and Wi, i = 1, . . . , M^, with Mk 
equal to the number of atoms in dictionary k, are normalization constants. Considerations 
are limited to the cases m = 2 (hat atoms) and m = 4 (atoms arising by discretizing cubic 
splines). Because we adopt the cut off approach for the boundary, the numbers Mk of 
total atoms in the fcth-dictionary varies according to the atom's support. For the linear 
spline based dictionaries the corresponding supports are: 1, 3, and 5, while for the cubic 
are 3, 7, and 11. 

With these dictionaries we construct the tensor product dictionary V = V^^Vj, i,j = i, . . . ,4 
which implies, approximately, a redundancy of five. However, the redundancy does not modify 
significantly the complexity order. The complexity of applying the OMP approach is dominated 
by the evaluation, at each iteration, of the inner products between the residual and the dictio- 
nary atoms. In the case of the proposed dictionaries the inner product with the DC dictionaries 
can be implemented by fast DCT and the complexity in computing the inner products with 
the other atoms depends on the atoms support (cf ([2])). By denoting di to the support of the 
spline based dictionary i, the complexity of computing the inner products at the selection step 
(0) is 

4 

0{{2Nf log, 2N) + 0{N J2 d^M^). 

i=2 

For the supported atoms we are considering the number of atoms Mj, i = 2,..., 4 is not 
much larger than A^. Moreover, for the hats dictionaries di = 1, = 3 = 5, so that the 
redundancy do not affect so much the stepwise complexity dominated by the inner products of 
the residual with all the dictionary's atoms. In addition which can be also accelerated by parallel 
calculations. Of course, the total complexity depends on the sparsity, as the complexity for 
selecting each atom has to be multiplied by the number of selected atoms, which is the feature 
of stepwise Pursuit Strategies. Nevertheless, the fact that the image is processed in small blocks 
leaves room for fast implementation by parallel processing. 

As can be observed in Table I, the performance in sparsity that is achieved with the pro- 
posed dictionaries is slightly better when using hats dictionaries, but with both dictionaries the 
sparsity in representing the six test images is about double than that yielded by faster nonlinear 
techniques such as DCT and WT. 



5 Conclusions 

Mixed DC and spline based dictionaries for sparse image representation have been introduced. 
It was shown that, comparing with fast nonlinear DCT and WT approaches, the proposed 
dictionaries yield a significant gain in sparsity. The complexity analysis and the fact that the 
processing is suitable for parallel computing lead to conclude that the proposed dictionaries can 
be of assistance to those image processing applications that benefit from the sparsity property 
of a representations. 

5 



Figure 3: The six test images from left to right, top to bottom: Boat, Bridge, Film chp, Lena, 
Mandril, Peppers 



Image 


DCT2 U Linear Splines 


DCT2 U Cubic Splines 


DCT 


Wavelets 


Boat 


7.05 


6.89 


3.63 


3.65 


Bridge 


4.24 


3.97 


2.06 


2.2 


Film 


9.72 


9.26 


4.53 


4.8 


Lena 


11.78 


11.7 


6.5 


6.97 


Mandril 


3.72 


3.5 


1.91 


1.90 


Peppers 


8.9 


8.62 


4.36 


3.39 



Table 1: Compression ratio (corresponding to PNSR=40 dB) achieved by each dictionary. The 
first column corresponds to the dictionary Vi composed of DC redundancy 2 and linear spline 
atoms of support 1, 3 and 5. The second column corresponds to dictionary V2 and cubic spline 
atoms of support 3, 7 and 11. The third column corresponds to the result obtained by nonlinear 
selection of DCT coefficients. The last column is the compression ratio produced by the Cohen- 
Daubechies-Feauveau 9/7 wavelet transform computed with the software WaveletCDF97 by 
thresholding coefficients so as to achieve the required PNSR or 40 dB) 
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